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A time-dependent Kohn-Sham (KS)-like equation for N bosons in a trap is gener- 
alized for the case of inelastic collisions. We derive adiabatic equations which are 
used to calculate the nonlinear dynamics of the Bose-Einstein condensate (BEC) 
and non-mean field corrections due to the three-body recombination. We find 
that the calculated corrections are about 13 times larger for 3D trapped dilute 
bose gases and about 7 times larger for ID trapped weakly interacting bose gases 
when compared with the corresponding corrections for the ground state energy 
and for the collective frequencies. 
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The newly created Bose-Einstein condensates (BEC) of weakly interacting alkali-metal 
atoms [1] stimulated a large number of theoretical investigations [2]. Most of these works 
are based on the assumption that the properties of BEC are well described by the Gross- 
Pitaevskii (GP) mean-field theory [3]. 

According to the theory of Lee, Huang and Yang (LHY) [4] the non-mean-field (NMF) 
correction to the GP theory behaves like y a 3 n, where n is the atomic density and a is the 
s-wave scattering length of interatomic interaction. The gas parameter, a 3 n, can have a large 
value when the number of atoms N in the BEC or a is large. In recent experiments [5], N 
was of the order 10 s with an intermediate value of the gas parameter, a 3 n « 10" 3 . 

Recently, it has become possible to tune the atomic scattering length to essentially any 
values, by exploiting the Feshbach resonances (FR) [6-8]. We note that the FR do not simply 
increase the gas parameter [9-13], and we do not consider the FR in this letter. 

Theoretical investigations of the NMF corrections to the ground state properties and to 
the collective frequencies of trapped BEC have already reported in the literature [14-18]. 

Inelastic collisions are an important issue in the physics of ultracold gases. The main 
goal of this letter is to consider a nonlinear dynamics of the BEC due to the three-body 
recombination. We calculate the corrections to the rate of the three-body recombination 
due to the NMF effects. These corrections are calculated in the large N limit and at zero 
temperature. 

Our starting point is the Kohn-Sham (KS)-like equation [19] for N interacting bosons in 
a trap potential V ext 

.^d^ h 2 „ 2t ta t d(ne(n)) 

where n(f,t) —\ *&(r,t) | 2 and e(n) is the ground state energy per particle of a uniform 
system. 

In Ref.[19], a rigorous proof is given to show that the KS-like equation correctly describes 
properties of the To nks- Girardeau gas in general time- dependent harmonic trap in the large 
N limit. 

The condition for the applicability of Eq.(l) demands that the function \1> changes slowly 
on distances of the order of the mean atomic separation n" 1 / 3 . Since the characteristic 
distance of changing ^ for the trapped gas is of the order of the radius of the condensate R, 
one gets the condition iV ss nR 3 3> 1. 

In order to take into account atoms lost by inelastic collisions (background collisions, 
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dipolar relaxation, three-body recombination, etc.) we model the loss by the rate equation 

= ~ J X(r,t)df, (2) 

where x(^) = J2i=i hn 1 gi(n) , n l gi is the local /-particle correlation function and k\ is the 
rate constant for the /-body atoms loss. For atoms in BEC, this rate constant is reduced by 
a factor of /!, which arises from the coherence properties of condensate [20, 21]. 
The generalization of Eq.(l) for the case of inelastic collisions reads 

..fiNf Ti 2 „ 2t Tr _ d(ne(n)) T .^v^, ; i / n T 

A lH = "2^ V * + V »'* + ^aP* " 'J (3) 

We will call Eq.(3) as dissipative KS-like equation. We note that Eq.(3) can be obtained 
from Eq.(l) using replacement of e(n) by e(n) — iT/2, where 



r/2 = ^2(fi/(2n))k l / x l - x gi {x)dx. 
i=i Jo 



It can be proved that every solution of the KS-like equation (1) is a stationary point of 
an action corresponding to the Lagrangian density 

ih d^/* d$ Ti 2 

where the asterisk denotes the complex conjugation. Indeed, the substitution of C , Eq.(4), 
into the following Euler-Lagrange (EL) equations 

SCp = OCq _ d dC _ dC = 
5^* ~ d^* dtd(dV*/dt) <9(W*) ~ ' 

(5) 

5C dC d dC _^ dC 



5ty dtd(dV/dt) 9(V*) 

gives Eq.(l) and its complex conjugate equation. 

Now for the dissipative KS-like equation (3) we write the corresponding Lagrangian C as 
a sum of two terms, a conservative one Co, Eq.(4), and nonconservative one CJ , C = Cq + C 
[22-24]. Now the EL equations read 

dC d dCp dCp S£f_ = Q 

dty* dtd(dV*/dt) d(W*) 5V* 

(6) 

<9£n d dC n _ dCa 5C 



(7) 



A comparison Eqs.(6) with Eq.(3) gives 

6¥~d¥~dtd(dV*/dt)~ <9(W*)~~2~^ lU ' 

sc _dc d dc dC _ih^ u 

dtd(d*/dt) y w^)~^h 9l{n) 

We are now ready to rewrite the Hamilton principle 5 J dt J df(C + C) = as 

5 J dtL + jJ2 k iJ dt J drn l - l gi (n){^*5^ - ^5^*) = 0, (8) 

where L — J C^dr. The variational formulation, Eq.(8), is an extension of the standard 
variational formulation for the situation where a Lagrangian corresponding to the original 
equations can not be found or does not exist [22-24]. 

For the remainder of this letter, we will focus and specialize solely on the three-body 
recombination. A suitable trial function can be taken as \I/(r, t) = exp[— (i/h)(j)(t)]Q(r, t), 
where both and <£> are real functions. With this ansatz, the Hamiltonian principle , Eq.(8), 
gives the following variational equations 

dN 



and 



a — = -^l 4> 6 (r, % 3 (4> 2 (r, t))df, (9) 



-<£ v2 * + + ^r 11 * = "«>*• (10) 



where N(t) = J $ 2 (f, t)dr, and fx(t) = d<p/dt. 

The condition of the validity of the adiabatic equations (9) and (10) is dN/dt < u> u N, 
where u u is a frequency of elementary excitation. 

The ground state energy per particle, e(n), in the low-density regime can be calculated 



using an expansion in power of Vna 3 

e(n ) = 2 J^ an [l + J^L( na 3)V2 + 8( ^ _ VS)na 3 [Hna 3 ) + C] + ...]. (11) 

The first term corresponding to the mean field prediction was first calculated by Lenz [25]. 
The coefficient of the (na 3 ) 3//2 term (the second term) was first calculated by Lee, Huang, and 
Yang [4], while the coefficient of the last term was first obtained by Wu [26]. The constant 
C after the logarithm term was considered in Ref. [27]. 
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The expansion (11) is asymptotic, and it was shown in Ref.[28] that the Lee- Huang- 
Yang (LHY) correction (second term in Eq.(ll)) represents a significant improvement on 
the mean field prediction up to a 3 n ~ 10~ 2 , but the logarithmic correction already is wrong 
at na 3 m 10~ 3 . In Refs. [16-18], the LHY expansion (first two terms in the expansion (11)) 
has been used to study effects beyond the mean field approximation. 

We do not consider the logarithmic term in the expansion (11), and rewrite Eq.(lO) in 
the limit of large N as 

m , _ x 32 



(H-V ext )=n{l + -^V^), (12) 

where n = $ 2 (f, t). 

At densities (na 3 < 1CT 3 ) Eq.(12) can be solved by iteration 

m . v 4m 3/2 , s/2 32am 2 . , 2 896a 2 m 5/2 . , 5/2 

(13) 

where 

„ km 1 / 2 ^? , 1024 9 ,a 2 ma TF ,25 22 , a 3 m 3 l 2 a%l , „, 

and the Thomas-Fermi (TF) approximation is simply obtained by keeping only the first term 
in the right side of Eq.(13), utf = mi^TF — V ext ) / (Anah 2 ) . 

Eq.(13) holds in the region where n is positive and n = outside this region. We note 
that the second term in Eq.(14) was considered in Refs. [14, 18]. 

For the harmonic trap potential V ext (f) = (m/2)(uj 2 x 2 + u 2 y 2 + uj 2 z 2 ), utf is given by 
[29] utf = (^ 2 /(2ma^ 5 )(15aiV) 2//5 , where a ho = (h/muj ho ) 1 ^ 2 is the oscillator length and 
uj ho = (oj x oj y uj z yi 3 . Eq.(13) in this case becomes 

(15) 



where utf(0) = itl^tf I '4:irah 2 . 

Using the correlation function [21] 



v 7F 

it can easily be seen that 

3357 



64 / 

ff3 (n) = (l + -=V^+4 (16) 



/x ■\-v\t 
&(r,t)g 3 (n)dr = -n 2 TF (0)N(l + — (Wn^O)) 1 / 2 + ...) (17) 
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Solution of Eqs.(9) and (10) in the mean-field approximation, corresponding to the first 
term in the expansion (16) and (17), reads 

t) - 8 ( a6/5a ^ /5 4. A h t)-vz m y n *\ 

$ M) " 21 { ^H0) + 5 ht) " A^tf VexU (18) 

and 

N(t) = (N-^(0) + l^y j tr 5/ \ (19) 
where a = 15 4 / 5 /1687r 2 , and for the rate r =| dY&N/dt | we obtain 

r = ^/a^alT (20) 

iV-4/5( ) + (4/5)(aA;3/a6/5af o /5 )t 

Analytical results, Eqs. (18-20), predict a strong Uho dependence of the 3-body recombina- 
tion in the TF regime, as shown for 87 Kb condensate with k% = 5.8 x 10~ 30 cm 6 / s [20] in Fig 1. 
For all cases of Fig.l, the adiabaticity is insured (N/N < 4x 10" 4 w feo and N/N < 5x lO^" 2 ^ 
for Who/2n = 12.83Hz and 77.78Hz respectivelly). 

A comparison of Eq.(17) with Eq.(15) shows that the non-mean-field corrections to the 
nonlinear dynamics of the BEC due to the three-body recombination are about 13 times 
larger than corresponding corrections to the ground state and to the collective frequencies 
of trapped BEC. As an example, we present in Fig. 2 the results of numerical calculation of 
the non- mean- field corrections to the rate r =| d\nN/dt | for the BEC of S7 Rb atoms with 
N(0) = 10 8 . 

For a ID Bose gas interacting via a repulsive 5-function potential, g5(x), the Lieb-Liniger 
(LL) model [30], e(n) is given by [30] e(n) = |^ 2 e(7), where 7 = mg/(fi 2 n) and for small 
values of 7, the following expression for e(n), e(n) = |(n — ^r^/^P + •••) is adequate up to 
approximately 7 = 2 [26]. In this case the large N limit solution of Eq.(10) reads 



n = 



(j,-V t 1 mg.fx-Vt 



+ -) 1/2 + -, (21) 



9 k\I h g 
where 

V = Vtf(1-\^) + ■■■), (22) 

with 7(0) = mg/h 2 n(0) and for V t = mu 2 x 2 /2, /i TF = (3Ngm 1/2 u/2 5/2 ) 2/3 . Using the ID 
correlation function [31] gz(n) = 1 — 6^/7/^ + we obtain 



24 973 1 

<S> 6 (x,t)g 3 (n)dx = -n 2 TF (0)N(l - —^7(0) + ...), (23) 



35 lpy ' v 512 
6 



where titf(Q) = Htf/q- 

The last equation (23) shows that for the weak interacting trapped ID bosons the NMF 
corrections to the three-body recombination dynamics are about 7 times larger than corre- 
sponding corrections to the ground state energy. 

In conclusion, we have developed a time-dependent dissipative Kohn-Sham (KS)-like 
equation for N bosons in a trap for the case of inelastic collisions. We derive adiabatic equa- 
tions which are used to calculate the nonlinear dynamics of the BEC due to the three-body 
recombination. The calculated non-mean-field corrections to the three-body recombination 
dynamics are shown to be about 13 times larger for 3-D trapped dilute bose gases and about 
7 times larger for ID trapped weak-interacting bose gases than corresponding corrections to 
the ground state energy and to the collective frequencies. 

We thank S. Khlebnikov for his interest and comments. 
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FIG. 1. The natural log of the number of 87 Rb atoms of the BEC in the trap as a 
function of time. The loss is due to the three-body recombination. Solid line and dashed line 
represent the theoretical results for traps with geometric average of the oscillator frequencies 
Uho/^TT = 12.83Hz and ushojlix = 77.78Hz, respectively. 
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FIG. 2. Rate of the three-body recombination, r =| |, as a function of time for the 
BEC of 87 Rb atoms for the initial condition N(t = 0) = 10 8 . The asymmetry parameter 
of the trap is A = v° and u z /(2ir) = 220Hz. The solid line corresponds to the mean- 
field approximation and the dashed line shows the results of inclusion the non mean-field 
corrections. 
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